Exact solutions for periodic and solitary matter waves in nonlinear lattices 
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We produce three vast classes of exact periodic and solitonic solutions to the one-dimensional 
Gross-Pitaevskii equation (GPE) with the pseudopotential in the form of a nonlinear lattice (NL), 
induced by a spatially periodic modulation of the local nonlinearity. It is well known that NLs in 
Bose-Einstein condensates (BECs) may be created by means of the Feshbach-resonance technique. 
The model may also include linear potentials with the same periodicity. The NL modulation func- 
tion, the linear potential (if any), and the corresponding exact solutions are expressed in terms of 
the Jacobi's elliptic functions of three types, cn, dn, and sn, which give rise to the three different 
classes of the solutions. The potentials and associated solutions are parameterized by two free con- 
' stants and an additional sign parameter in the absence of the linear potential. In the presence of 

' the latter, the solution families feature two additional free parameters. The families include both 

sign-constant and sign-changing NLs. Density maxima of the solutions may coincide with either 
minima or maxima of the periodic pseudopotential. The solutions reduce to solitons in the limit of 

Qthe infinite period. The stability of the solutions is tested via systematic direct simulations of the 
GPE. As a result, stability regions are identified for the periodic solutions and solitons. The periodic 
patterns of cn type, and the respective limit-form solutions in the form of bright solitons, may be 
stable both in the absence and presence of the linear potential. On the contrary, the stability of the 
two other solution classes, of the dn and sn types, is only possible with the linear potential. 
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PACS numbers: 



I. INTRODUCTION 



Due to the well-known fact that the Gross-Pitaevskii equation (GPE) provides for an exceptionally accurate descrip- 
tion of the dynamics of Bose-Einstein condensates (BECs) in rarefied ultracold gases |5fl] , and the related nonlinear 
J> ' Schrodinger equation (NLSE) is an equally good model of the light transmission in nonlinear media [s^ , finding exact 
, solutions to these equations in nontrivial settings is a problem of great significance [l^ [6l| . A large class of exact 
periodic solutions, which provide insight into the structure of one-dimensional (ID) BEG patterns, was found in the 
framework of the GPE with specially designed lattice potentials [lCll4l2l llil [TEl ISOl ISq , that can be created by means 
of available experimental techniques (2^ . Exact solutions were also found for pro pag ating density waves [25j , various 
localized states supported by particular spatiotemporal settings 0, H, [l^, [l9l l2ll. [371 138l.l40l |48| , |63| , stationary local- 
' ized modes supported by the self-attractive nonlinearity represented by sets of one or two delta- functions [4l|, [g^] , 
' and for settings providing for the reduction of 3D states to the ID form [23j . 
. . , Another noteworthy class of exact solutions was recently reported in a model combining linear and nonlinear lattices 
. ^ ' (NLs) In the experiment, linear periodic potentials can be induced by means of the well-known optical-lattice 
k>( ' techniques 0, [13, [s^ i46|] , while NLs represent spatially periodic modulations of the nonlinearity coefficient, which 
j_j ■ may be implemented by means of accor ding ly structured external fields via the Feshbach-resonance effect. In fact, 
Ci ■ exact solutions reported in Refs. [4l|, [4^ |62| and [1] are particular representatives of a vast class of states supported 
by various patterns of the spatial modulation of the nonlinearity (alias nonlinear pseudopotentials 1311 . in terms of the 
condensed-matter theory), which were investigated in ID |l[, ll El, [M [M \M, [5l[Hj58laiid 2D [H, [H [M IM Hi 
settings; mathematical aspects of these solutions were also studied in some detail (711281 l49l [57|. 

An issue of obvious interest is to find classes of exact periodic solutions supported by periodic pseudopotentials, and 
also ultimate cases of such states with the infinite period, representing matter- wave solitons. In addition, a challenging 
problem is to find exact solutions that may be supported by combinations of a nonlinear periodic pseudopotential and 
its linear counterpart of the optical-lattice type. Of course, these states may be physically meaningful only if they 
are stable. The objectives of the present work is to report three classes of exact periodic solutions for the GPE with 
a periodic pseudopotential, as well as for a more general model including a linear periodic (lattice) potential, and to 
analyze the stability of the periodic and solitonic solutions. We find the form of the periodic functions accounting for 
the shape of the NL pseudopotential and linear potential (if the latter is included in the model), and, simultaneously, 
the respective exact solutions by means of methods based on the Hirota bilinear representation, which were elaborated, 
in the context of "management" models for solitons, in Refs. [l9[, [2l|, [4^, [3, and for periodic (cnoidal) waves in 
various models with the bimodal (two-component) structure and complex nonlinearities in Refs. [lo, [13, [SO] ■ Making 
use of the basic Jacobi's elliptic functions, cn, sn, and dn, we are able to find three distinct families of nonsingular 
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solutions, affiliated with each of these functions; accordingly, the families are referred to as ones of the "cn" , "sn" , 
and "dn" types. In the case when only the NL pseudopotential is included, while the linear potential is absent, each 
family is controlled by three free parameters, two continuous and one which determines the sign of the nonlinearity. 
If the linear-lattice potential is included too, each family features four free continuous parameters. In the limit of 
the infinite period, the solutions of the cn and dn types go over into bright solitons, while the respective limit of the 
sn-type solutions is a dark soliton. The stability of all these solutions, both periodic and solitary (except for the dark 
solitons), is tested in this work by means of systematic direct simulations of the perturbed evolution. 

A noteworthy feature of the exact solutions of all the three types is that they can be found in cases when the local 
nonlinearity periodically changes its sign, thus forming NLs with the spatially periodic alternation of self-focusing 
and defocusing layers. It is relevant to note that exact solutions which were recently reported in a model combining 
nonlinear and linear periodic potentials 4] cannot be obtained in such they are generated by a gauge 

transformation of the NLSE with constant coefficients. 

The normalized form and physical meaning of the GPE for the atomic wave function, ^'(a;, t), and of the NLSE for 



the local amplitude of the electromagnetic wave in optics, are well known [36|,|5Qj (in the latter case, t is t he prop agation 
distance of the beam launched into the layered medium in the longitudinal direction, see Refs. [sl. ITtI. [26L H. [59l| ) : 

+ (1/2)^',, + g{x) * - y(x)* = 0. (1) 

Here g{x) and V{x) represent the nonlinear and linear lattices, respectively, with g{x) > 0/g{x) < corresponding 
to the local attraction/repulsion between atoms in the BEG. As usual, stationary solutions to Eq. ([Ij with chemical 
potential are sought for as ^' (x, t) — ip{x) exp(— z/ii), where real function ip{x) obeys equation 

fi^P+{l/2)ij" + g{x)ij^ -V{x)ij = 0. (2) 

The corresponding expression for the energy of the matter-wave configuration is 



E 



1 ii/f - ^gix)^^ + Vix)^' 



dx. (3) 



The analysis presented in this work is focused, first, on looking for special forms of periodic functions g{x) and V{x), 
for which Eq. ^ admits exact solutions in terms of elliptic functions. Then, the stability of the exact solutions is 
tested by means of simulations of the perturbed evolution of those solutions within the framework of Eq. 

The paper is organized as follows. In Section 2, we introduce the class of exact solution of the cn type, first in the 
model without the linear potential, and then in the its more general version, including the potential. Results of tests 
of the stability of this class of exact periodic solutions are reported in Section 3. The cn patterns may be stable both 
without the linear potential, and in the presence of the potential. Solutions of the two other classes, of the dn and 
sn types, together with the results of tests of their stability, are presented in Sections 4 and 5, respectively. Both 
these species of the patterns require the presence of the linear potential for their stability. Section 6 reports results 
for bright solitons, which can be obtained as ultimate-form solutions from families of the cn and dn types. The bright 
solitons may be stable in the absence and presence of the linear potential alike. Dark solitons, which are limit-form 
solutions of the sn type, are briefly considered too in Section 6, without the stability analysis, which will be reported 
elsewhere. Conclusions drawn from this work are formulated in Section 7. 



II. EXACT SOLUTIONS OF THE CN TYPE 

A. The model without the linear potential 

We first consider the case of the NL pseudopotential alone, with V{x) = in Eq. ([1]). Using the techniques 
elaborated in Refs. [60j . which are based on manipulations with the Jacobi's elliptic functions, it is 

possible to identify general forms of the NL structural function, i.e., coefficient g{x) in Eq. ([J), which admit exact 
solutions to Eq. ^ with V = 0, and the form of the solution, ip{x), itself. The first type of the solutions is based on 
the elliptic cosine (cn), the corresponding expressions for g(x) and tpix) being 

- ':ir^^'t (4) 

1 + cn^(ra;) 

V'(x)=Ao— = ; ; (5) 
y'l + o cn"^(rx) 
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The spatial periods of solution ([5]) and NL potential ([4]) are, respectively, 

L = AK{k)/r (6) 

and where k is the modulus of cn, which takes values < fc < 1, K{k) is the complete elliptic integral of the first 
kind (the case of fc = 1 corresponds to L = oo, i.e., localized solutions). Generally, expressions Q and ([5]) contain a 
set of four free parameters, viz.^ go-, r, fc, b, while the remaining coefficient (modulation depth), gi, together with the 
amplitude and chemical potential of the solution, are found to be 

_ ffob b(l + 3&) + (l + b)(l-3b)fc^ 

2 6(2 + 36)- (l + fe)(l + 36)fc2' ^ > 



= .gor^[-6(2 + 36) + (1 + 6)(1 + 36)fc^], (8) 

(rV2) [1 + 36- (2 + 36) fc2]. (9) 

Using the obvious scaling invariance and assuming go ^ 0, we fix |go| = 1 (if go — 0, the above solution becomes 
trivial). Throughout this paper, (70 = il is kept as the sign parameter. Thus, fc and r determine the period of the 
NL (L), go is the overall sign of the nonlinearity, and 6 controls the depth of the spatial modulation of the periodic 
pseudopotential. An obvious constraint on 6, which may be both positive and negative, is 6 > — 1, which is necessary 
to avoid singularities in Eqs. Q and ([5]). Further, the remaining scaling invariance of Eqs. ([Ij and ([2]) makes it 
possible to set r = 1, without the loss of generality, which is fixed below, unless the linear potential is included. Thus, 
in the absence of the linear potential, the family of the exact solutions of the cn type depends on two free coefhcients, 
fc and 6, and the sign parameter, go. 

Comparing expressions Q and ([5]), it is easy to see that maxima of the density of the periodic solutions, ip'^{x), 
which are always collocated with points where cn^irx) = 1, coincide with maxima or minima of g{x)-i.e., with, 
respectively, minima or maxima of the NL pseudopotential-in cases when, severally, gi > gob or gi < gob. Then, 
substituting expression ([7]) for gi, we conclude that, for 6 < 0, the density maxima always coincide with maxima of 
g{x), while for 6 > this is true under condition 

fc2<6/(l + 6). (10) 

Otherwise, density minima {ip = 0) are located at maxima of g(x). In fact, taking into account condition Aq > 0, as 
it follows from Eq. ([8]), one can conclude that condition (fTO)) may only hold in the case of go — ^1- Further, it is seen 
from Eq. (jS]) that the ground state, which realizes the minimum of the energy, must have density maxima collocated 
with minima of the pseudopotential, i.e., maxima of g{x) (of course, this condition is only necessary but not sufficient 
for the exact solution to represent the ground state). It is relevant to mention that various localized modes in the 
simplest model of the NL, represented just by two spots at which the self-attractive nonlinearity is concentrated, may 
be stable, even if they do not represent the ground state f43|. 

Lastly, we note that, as it follows from expression Q, Eq. ([2]) degenerates into the ordinary GPE, with g{x) = ±1, 
in two cases: 6 = 0, or 

91 = gob. (11) 
In the former case, expressions ([S]) and ^ go over into the usual cnoidal solution of the GPE with g{x) = +1, i.e., 

^-(2;) = fc cn(a;), ^ = (1/2) (l - 2k^) , (12) 

provided that 50 — +1 [for go ^ —I and 6 = 0, the solution does not exist, as Eq. ^ yields A§ — —k^]. The other 
degeneration condition, given by Eq. pT|) . if combined with Eq. (O, yields 

fc2=6/(l + 6). (13) 

In this case, solution ^ is relevant for go = —1, being 

tjj(x) = fc cn(a;)/dn(a;), (14) 

[for go = +1, one obtains Ag < from Eq. (|8|)]. It is easy to check that expression ([T4| is a straightforward exact 
solution to Eq. ^ with g{x) = —1 and ^ = (1/2) (l + fc^), in agreement with Eq. (0). 
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B. The model with the linear potential 



A more general class of exact periodic solutions of the cn type can be found if the NL pscudopotential, introduced, 
as above, in the form of Eqs. (U) [but Eq. ([7]) should be dropped, see an explanation below], is combined with the 
following linear-lattice potential: 

_ l + yocn^(rx) 

^^''^^ l + 6cn2(ra;)' ^^^^ 



Vo = b + —[b{l + 3b) + (1 - 3b){l + b)k^] 

+ |^^^[6-(l + 6).^]. (16) 
2(51 - gob) 

The respective exact cnoidal solutions keeps the same general form as in Eq. ([S]), but with the amplitude and chemical 
potential given by 

3br^(l + 6)[b-(l + 6)P] 
2 [91 - gab) 



+ (rV2)[(l + 36) - (2 + 3b)k^]. (18) 

This solution family depends on four continuous real parameters, viz.^ gi,r,k,b^ and sign parameter go. Unlike the 
case of y = 0, coefficient r cannot be scaled out, unless 6 = 0, and the value of gi is now another free parameter, 
rather than the one given by expression Q; in fact, Eq. ([7]) follows from Eq. (fT6|) if one sets Vb = 6, which reduces 
linear potential (|15l) to a trivial form, V{x) = 1. 

In the limit of 6 — > 0, one obtains, from the above formulas, g{x) = go and V{x) = 1 + (1/2) (rfc)^ cn^(rx), the 
respective exact solution ([5]) being ip{x) = Aocn{rx) with An = ( 3/2)gn ( rk) ^ and fi = 1 — (fc^ — 1/2), provided 
that go = +1. This is one of exact solutions found in Refs. [lO, [U, 13 US HP] for Eq. ^ with constant g. 



C. Further analysis of the exact solutions (without the linear potential) 

The above solutions are meaningful provided that Eqs. ([8]) and pTl) yield Aq > 0. Because the exact solution is 
defined in a vast parameter space, several cases should be considered separately. In this subsection, we do that for 
the solutions obtained in the absence of the linear potential. 

The first particular case corresponds to go = +1, 6 > 0, see Eq. (j?]). Then, as follows from Eq. ([5]), condition 
Aq > requires 

k^ > ^(^ + ^^) (19) 
> (l + 6)(l + 36)- ^'^^ 



For all 6 > 0, Eq. (|19p is compatible with constraint fc < 1. On the other hand, Eq. p9|) is incompatible with 
condition (jlO[) . hence this subfamily of the exact solutions cannot represent the ground state (and, in fact, it cannot 
be stable, as shown below). 

As said above, an interesting possibility is to find exact solutions in the sign-changing NL, which makes the model 
drastically different from the usual NLSE. As follows from Eqs. @ and ([7]), the NL-modulation coefficient, g{x), 
periodically changes its sign for gi < —1, i.e.. 



6(l + 36) + (l + 6)(l-36)A:2 



6(2 + 36)-(l + 6)(l + 36)fc2 



^ < -1- (20) 



In view of underlying condition ([T9|) . the denominator on the left-hand side of Eq. pO|) is negative, hence this 
inequality, together with Eq. (|19p . identifies a parameter region in which the exact solutions are available for the 
sign-changing NL: 

6(2 + 36) K4 + 36) _ ^^^^ 



(l + 6)(H-36) (l + 6)(2 + 36)' 
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It is easy to check that, for 6 > 0, interval (|2T|) always exists, and, as a whole, it lies within the region of < fc < 1. 

At small values of 6, region (PT|) is squeezed into a narrow strip, 2b— 56^ < < 2b— 76^/2, where NL-modulation 
function ([4]) takes the form of g{x) ~ 1 + gi cos^ x [recall that, while considering Eq. ([1]) with = 0, we set r = 1 
in all expressions following from Eqs. ^ and ([5])], and exact solutions ([5]), (|S]) become quasi- harmonic ones, with a 
small amplitude, ip{x) « by^—3/ {2gi) cosx, the corresponding chemical potential being n « 1/2. In this asymptotic 
form of the solution, gi may be treated as an arbitrary parameter ~ 0{1), the solution existing solely for gi < 0, 
which includes the case of the sign-changing NL for (/i < — 1. It is evident that this explicit asymptotic solution has 
its density maxima coinciding with minima of g{x) (maxima of the pseudopotential) , which are located at a: = nn, 
n = 0,±1,±2,... . 

In the opposite limit, 5^1, existence region (PT|) is narrow too, 1 — 2/ (36) < fc^ < 1 — 1/ (36). The respective 
modulation profile, g{x), and exact solution ([5]) reduce to chains of "flat-top" solitons. 

The case of go = +1 and — 1 < 6 < in the nonlinearity-modulation function (|4]) can be considered too (recall the 
solutions with 6 < may represent the ground state). The inspection of Eq. ([8]) shows that exact solution ([5]) exists 
for ah fc2 < 1 at < -6 < 1/2, and, at 1/2 < -6 < 2/3, it exists in interval 

, |6|(2-3|6|) 
<(3|6|-1)(1-|6|)' (''^ 

cf. Eq. p^ . there being no solutions at — 1 < 6 < —2/3. Further considerations of Eqs. (|3]) and (0 demonstrate 
that the NL cannot change its sign in this case, always corresponding to the local self-attraction. 

The other generic subfamily of the exact solutions is obtained for go = -1- First, considering this case with 6 > 
in Eq. Q, condition Aq > imposes a constraint on the elliptic modulus, 

6(2 + 36) 

< (l + 6)(l + 36)' ^^^^ 

which always complies with fc < 1. Further straightforward considerations demonstrate that constraint ([25)) does not 
admit sign-changing NLs, i.e., in the case of go = —1 and 6 > 0, the exact solutions pertain to the spatially modulated 
nonlinearity which remains self-repulsive at all x. 

Another possibility is to take 5o = " 1 a-nd — 1 < 6 < in Eq. (U). A simple algebra demonstrates that the 
corresponding condition Aq > [see Eq. ([5])], combined with fc < 1, is satisfied with any k for 2/3 < —6 < 1, and is 
never satisfied for < — 6 < 1/2. In the intermediate case of 1/2 < —6 < 2/3, the exact solutions exist in interval 

H (2-3161) 

>(3|6|-1)(1-|6|)' ^2'^ 

cf. Eq. (|22l) . The consideration of expression ([7]) demonstrates that condition gi > 1 automatically holds in the 
present case, hence NL modulation function @ is always of the sign-changing type. Thus the exact solutions in this 
case may, in principle, represent the ground state of the BEC with nonzero mean density, loaded into the sign-changing 
NL. In fact, it will be demonstrated below that solely in this case the stationary periodic solutions of the cn type may 
be stable (without the addition of the linear potential). 

It is relevant to stress that the existence of the exact cn-type solutions with 170 = — 1 is a nontrivial finding in the 
following sense: according to Eqs. (O, (O, (|S]), and ©, in the limit of 6 = the same type of the solution, but 
taken with go = +1, goes over into the commonly known exact unstable cnoidal solution for the GPE with g{x) = +1 
(the uniform self- attraction), given by expressions (jl2p . Thus, the cn solutions for go = +1 may be regarded as a 
continuation of this simple unstable solution, which helps to explain the fact that they are always unstable too, as 
shown below. On the other hand, in the case of go = —1 the existence range of the exact solutions is separated, 
as demonstrated above, by interval < — 6 < 1/2 from 6 = [in particular, at 6 = and go = —1, Eq. dH) yields 
Aq = — fc^], i.e., this subfamily of the cn solutions has no counterpart among solutions of the GPE with the constant 
nonlinearity coefficient, and it does not originate from that limit-in particular, it may be stable for this reason. 



III. THE STABILITY OF THE CN-TYPE SOLUTIONS 



The test of stability of the exact solutions is a crucially important part of the analysis. In this work, we do 
not aim to compute stability eigenvalues for small perturbations around stationary solutions. Instead, we focus on 
systematic direct simulations of perturbed solutions. Although less rigorous mathematically, this approach is closer 
to the description of the physical situations, where finite random perturbations are always present. 
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FIG. 1: (a) The spatiotemporal evolution of the density of a numerical solution to Eq. ([l]) demonstrates the onset of the 
instability of the exact periodic solution of the cn type, induced by a weak initial random perturbation (at the level of 1% of 
the amplitude of the stationary pattern), (b) The corresponding profile of the periodic nonlinearity-modulation function g{x). 
The respective values of parameters in Eqs. Q, ([5]), and (|6]) are go = +1, fe = 1, r = 1, and k = 0.82. This exact solution was 
obtained in the absence of the linear potential, with gi — —2.25 given by Eq. 0. 

The numerical integration was performed in the spatial domain with the length corresponding to 100 periods of 
the structure, see Eq. and periodic boundary conditions. The simulations were run up to time t = 5000, in 
terms of the present notation. For most patterns, which have the spatial period L < 5 (see below), this implies 
t > 100 characteristic dispersion times, estimated as Tdisp ~ L^- In fact, in all cases the instability, if any, manifests 
itself at a much shorter time scale, t < 500. The simulations were extended to t = 5000 to make sure that solutions 
identified as a stable ones do not develop a delayed instability (still longer simulations, carried out for selected cases, 
completely corroborate the stability). Initial random perturbations, at the 1% amplitude level, were always sufficient 
to unambiguously categorize stable and unstable solutions (the application of stronger perturbations did not change 
the results). 

A. The solutions without the linear potential 

We start the presentation of the stability results by considering the cn-type solutions in the absence of the linear 
potential, V{x) = in Eqs. ([IJ and ([2]). First of all, in all cases with go = +1 and b > (which suggests that 
the nonlinearity is self- attractive, on the average), the simulations reveal that the cn-type solutions are unstable, see 
a typical example in Fig. [TJ In accordance with results of the above analysis, density maxima of the unperturbed 
solution in Fig. [TJa) coincide with minima of g{x) in Fig. [Ijb), i.e., maxima of the respective pseudopotential, which 
readily explains the instability. Note that, although only six spatial periods of the stationary solution are shown in 
the figure, it is actually a representative fragment of the picture which was generated in the domain covering 100 
periods, as said above. The same pertains to other figures displayed below. 

A generic feature observed in Fig. [T]is that the instability sets in after evolution time t < lOTdisp. The instability 
onset may be delayed for the unperturbed solutions with a small amplitude, which is an obvious manifestation of the 
nonlinear nature of the instability. 

In the same case of qq = +1 but with 6 < 0, unambiguously stable stationary solutions were not found cither. 
The difference of this case from the one with & > is that the respective NL function g{x) is a sign-constant one, 
i.e., the nonlinearity is locally attractive everywhere (which makes the modulational instability of stationary periodic 
solutions quite plausible). 

Proceeding to the case of go = ~lj i-e., the model with the nonlinearity which tends to be self- repulsive on the 
average, no stable solutions have been found for 6 > in Eq. (jH). However, stable solutions are possible if 50 = —1 is 
combined with & < and appropriate values of fc, which should be relatively close to fc = 1 [recall that the respective 
exact solutions exist in the interval of 1/2 < —6 < —1, and they all correspond to a sign-changing function g{x)]. 
An example of a stable pattern is shown in Fig. [21 In fact, this figure displays a marginally stable solution, in the 
sense that it is never destroyed by the originally imposed random perturbation, but the disturbance does not dissipate 
either, remaining trapped in the solution. Note that this example shows density maxima collocated with maxima 
of g{x), i.e., minima of the corresponding NL pseudopotential, which facilitates the stabilization of the pattern, and 
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(a) (b) 

FIG. 3: The same as in Fig. O but for k = 0.70 and the corresponding value of gi = 1.85, as given by Eq. ([7]). In this case, 
the cn-type solution is unstable. 

gives it a chance to be a ground state. 

Nevertheless, in the case oi go ~ ~1 and 6 < the cn solutions may be unstable if the elliptic modulus is taken too 
far from k = 1. For instance. Fig. [3] displays an example of a quick onset of the instability, found at the same value 
of 5 = -0.7 as in Fig. [21 but with k = 0.95 replaced by fc = 0.70. 

Collecting data of the systematic numerical simulations for the presently considered subfamily of the exact cn-type 
solutions, with go = —1, 1/2 < —b < 1 and various values of fc, makes it possible to identify a stability area in the 
plane of (6, fc), as shown in Fig. 21 The stability boundary in this figure was drawn by interpolating results obtained 
from the direct simulations for 

&=-0.9, - 0.85, - 0.8, - 0.75, - 0.7, - 0.67, - 0.65, - 0.6, - 0.55, - 0.51. 

As said above, the stability is limited to relatively long-wave patterns, with k sufficiently close to 1, viz., k > 0.9, 
the respective values of period ([6]) being L > 10.12. In accordance with the above analysis, the existence region of 
the solutions is limited to 1/2 < —6 < 1. It may be relevant to notice that the stability boundary in Fig. [His not 
monotonous, which may be understood as a consequence of the fact that the solution depends on 6 in a complex 
manner, as seen from Eqs. (O and 
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FIG. 4: The region of the stabihty of the cn-type periodic stationary solutions, as found for go — —1 without the hnear 
potential, in the plane of parameters b and k, which are defined in Eq. Q. The dashed oblique line with arrows indicates the 
border of the area where the solutions exist, as per Eq. (|24p . for 1/2 < —b < 2/3. 

B. The cn-type solutions in the presence of the linear potential 

The linear potential given by expressions ((T5]) and affects not only the shape of the solutions but also their 
stability. First of all, in the cases of go = +1 with cither sign of b, and go — —I with 6 > 0, well-defined stable 
solutions have not been found, as well as in the same cases in the absence of the linear potential, see above ("well- 
defined" implies that we do not count solutions with a very small amplitude, that may seem stable as the respective 
nonlinearity is extremely weak) . A typical example of a retarded onset of the instability of a small-amplitude solution 
is displayed in Fig. [Sj for go = +1 and b = +0.5. 

Similar to the situation reported above for the solutions with V{x) = 0, a vast subfamily of stable solutions is 
found in the case of go ~ —1 and 6 < 0. Actually, the stability region is located at 6 > —0.73, this border being 
practically independent of elliptical modulus fc (at 6 = —0.72, the exact solutions are found to be stable in interval 
0.01 < k < 0.99, while at 6 = —0.74 they are unstable in the same interval). The shape of the stable solutions is 
quasi- harmonic at small and moderate values of k (see an example in Fig. |6]), while featuring sharp peaks at k close 
to 1 (Fig. [7])- Note that the amplitude of the solutions shown in these typical examples is not especially small, which 
confirms that the observed stability is not a trivial consequence of the weak nonlinearity. It is also worthy to note 
that density maxima of the solutions coincide with maxima of the NL modulation function, 5 (a:), and minima of 
y(a;), i.e., the density maxima are collocated with minima of both the nonlinear pseudopotential and linear potential. 
This circumstance obviously facilitates the stabilization of the periodic patterns, giving them a chance to realize the 
ground state. 

IV. EXACT SOLUTIONS OF DN TYPE AND THEIR STABILITY 
A. Stationary solutions 

Another generic class of exact solutions, which, as well as the solutions of the cn type, are represented by even 
functions of a;, can be obtained from cn by the continuation to values of fc > 1, using the well-known transformation 
formula, 

cn{x,k) = dn{kx,l/k) . (25) 

Thus, in the absence of the linear potential, V{x) — 0, expressions generate the following class of exact solutions 

to Eq. ©: 




10 
(a) 




(b) (c) 



FIG. 5: Panels (a) and (b) have the same meaning as in Figs. [T]-[3l but here we display an example of the perturbed evolution 
of a weakly unstable cn-type solution, found in the presence of the linear periodic potential, as given by Eqs. (|15p and (|16|) . 
The shape of the linear potential is displayed in panel (c). The parameters are go = +1, b — +0.5, r — 1, gi — —2, and 
k = 0.60. 

^{x) - A. ^ ' =, (27) 

Jl + b dn^irx) 



^ 9ob {b + l){3b-l)-b{3b+l)e 

2 (fe + 1) (36+1)- 6(35 + 2)P' ^ ^ 

Al = .gor^ [(5 + 1) (36 + 1) - 6(36 + 2)k^] , (29) 

/i= (rV2)[-2 - 36+(36+l)fc2], (30) 

where, as above, we assume normalization go — ±1, the elliptic modulus takes values < k < 1, and an additional 
normalization may be imposed, r = 1. 

As concerns the transition to the ordinary GPE, with g = const, Eq. (f26| demonstrates that this is possible in 
two cases, 6 = or under condition coinciding with Eq. PT|) . see above. The latter condition, if applied to Eq. (|28p . 
leads to fc^ = (6 + 1) /6. Unlike the similar condition obtained for the cn-type solutions in the form of Eq. ([T3|. the 
present one, obviously, cannot give appropriate values of the elliptical modulus, < fc^ < 1, if 6 takes values for which 
solution (l?71) is nonsingular, i.e., 6 > — 1. As for the former limit case, 6 = 0, expressions (P7)) and yield the 
ordinary unstable solution to the equation with the uniform self-attraction, g{x) = +1, in the form of ip{x) = dn(a;). 

Making use of Eqs. (|T5|) - (|T8l) taken at fc > 1 and reducing this to values fc < 1 by means of Eq. (|25l) . the dn solution 
family may be extended to the case when the linear potential is present in Eq. 



(31) 
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j,2 Ijj,2 



Vo^-+"—[{3b+l)e-{3b + 2)] 



3(6+1)51 ^^,2 



[bk'-{l + b)], (32) 



_3br^il + b)[bk^-{l + b)] 

2 {gi - gob) 



^=(rV2)[(36+l)A:2-(2 + 36)]. (34) 

Similar to the class of the en-based solutions, in the presence of the linear potential the NL-modulation function is 
given by the same general expression as above, Eq. ([26l) . while Eq. ([28|) should be dropped (actually, the latter 
equation is the condition for reducing Vq to zero). Accordingly, the solution family depends on four free parameters, 
gi,r, k, b, and the sign parameter, go- 



B. The stability of the dn-type solutions 

In the absence of the linear potential, the exact solutions given by Eqs. are found to be strongly unstable, 

see a typical example in Fig. [8] Also unstable are the exact solutions found in the presence of the linear potential, 
as per Eqs. (|3T |) - (p3l) . in the case of 50 = on the contrary to the solutions of the cn type, and also to the odd 
solutions of the sn type (see below), where the potential can readily stabilize the solutions with go = —1. 

On the other hand, in the case of go = +1, when the solutions of both the cn (see above) and sn (see below) types 
cannot be stable, the dn species of the periodic patterns readily demonstrates its stability, in the presence of the 
linear potential, as shown in Fig. [S] Note that panels (b) and (c) of this figure imply the competition between the 
nonlinear pseudopotential and linear potential, as minima of g{x) coincide with minima of V{x), this sets of points 
also coinciding with density maxima of the solution, see panel (a) . A similar competition between the nonlinear and 
linear effective potentials was recently considered in a model of a one-dimensional photonic crystal [l^l- In fact, almost 
all the solutions belonging to this subfamily are stable, with the exception of ones with large 6 or fc very close to 1. 



V. SOLUTIONS OF THE SN TYPE AND THEIR STABILITY 



A. Solutions without the linear potential 



In addition to the two classes of even stationary solutions considered above, it is possible to find a class of the 
NL-modulation functions, g{x), and related potentials, V{x), which are even too, but they support exact stationary 
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FIG. 9: The same as in Fig. [S] but for a typical stable solution of the dn type found in the presence of the linear potential, for 
go — +1, b — 1, r — 1, gi — —2, and k = 0.9. 



solutions which are odd in x, that are based on the Jacobi's function sn. First, in the case of V{x) ~ 0, the explicit 
form of these solutions, along with function g{x) which supports them, is [cf. similar expressions (jll)-® and p6|) -(|30 |) 
for the solutions of the cn and dn types] : 



go +gisn^(ra;) 
1 + 5 ird^irx) 



(35) 



74osn(rx) 
i/l + 6 sn^(ra;) 



(36) 



go6 6(1 + 36) - (1 - 6)fc^ 
~6(2 + 36) + (l + 26)fc2 



(37) 



Al = -gor2[6(2 + 36) + (1 + 2h)k\ 



(38) 



^ = _ (1 + 36 + A:^) 



(39) 



In the case of go = +1, expression (p8| cannot give > ^ iox h > Further, a straightforward analysis 
demonstrates that, in the same case but for 6 taken from < —6 < 1/3, the solutions exist in interval 



(40) 
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of values of the elliptic modulus. This entire interval satisfies constraint fc < 1. For 1/3 < —b < 2/3, the solutions 
exist for all < 1. Finally, for2/3<— 6<1, the existence condition is 

IH(3|6|-2) 



2|6|-1 



< fc^ < 1. (41) 



A similar analysis performed for the case of go = —1 demonstrates that, with & > 0, the solutions exist for all values 
of fc^ < 1. Proceeding to negative b in this case, it is easy to conclude that, for 1/3 < — 6 < 0, the solutions with 
> exist for 

l&l(2-3|6|) ^ 



1-2|6| 



< fc" < 1, (42) 



cf. existence condition (j40p for go = +1, which was obtained in the same interval of b. In the adjacent interval, 
1/3 < — 6 < 2/3, the solutions do not exist at all. Finally, in the end portion of the allowed interval of b, viz., 
2/3 < —b < 1, the sn-type solutions are available in the region of 



,2 . \b\ m - 2) 



' < (43) 

cf. region (PT|) in the case of go = +1 [the whole interval (|33]) complies with constraint fc < 1]. 

The transition of g{x) given by Eq. ([35]) to g = const, corresponding to the ordinary GPE, may be performed 
by taking & = 0, or by imposing the same condition (jll[) as above. In the former case, the respective solution, 
ip{x) — k sn(a;), satisfies the GPE with g{x) = —1. In the latter case, the solution is relevant for g{x) = +1, reducing 
to 6 = — fc^, ^ = (1/2) (l — 2fc2), and i^{x) — k^/1 — k'^sn{x)/dn{x), cf. similar hmit solution (fM)) of the cn type. 

B. Exact sn solutions with a linear potential 

Another family of the sn solutions can be found, adding to Eq. ([2]) the linear potential in the following form, 

1 + sn^(ra;) 

fcV^ 6(l + 36 + fc2)r2 3gi{l + b){b + k^y , , 

""'^-2 2 2ig,-gob) ' (45) 

Then, ip{x) is again given Eq. ([36)) . but Eqs. (l38l) and ([39]) are replaced by 

2 _ 3b{i + b){b + k^y 

2{gi - gob) 



/i= (rV2)(l + 3fe + fc2), (47) 

while Eq. ([57]) is dropped. 

In the limit of 6 — > and — > 0, the present solution makes sense for go — —1- In this limit, we may again fix 
r = 1, hence Eqs. ([M]), (O and ^ yield 

1 + fc2 

ip = V3/2fc sn(x), ^l = (48) 

Expressions (|48|) provide for a solution to an equation with the constant nonlinearity coefficient, 

fxtp + (1/2)1/;" -tP^ + (fcV2) sn^{x)i! = 0. (49) 
These particular results, obtained for g(x) = —1, reproduce findings reported in Ref. [12i] . 
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C. The stability 

In the absence of the hnear potential, nearly all the solutions of the sn type are unstable, as illustrated by a typical 
example in Fig. (TU] Note that, in this case, maxima of the density coincide with minima of g{x), thus pushing the 
solutions towards instability. The only stable solutions could be found for go — —1 and very small when, as said 
above, the exact solutions are close to the commonly known stable exact solution, = k sn(a;), to the GPE with 
9{x) = -l. 

While, unlike the cn solutions, all patterns of the sn type tend to be unstable without the linear potential, in the 
presence of the potential they may be readily stabilized, in the case of go = ^1 smd b < 0, which resembles the 
respective property of the cn-type patterns, see Fig. U) An example of the stable sn structure is shown, for this case, 
in Fig. [TTJ As is typical for other types of stable patterns, in this figure we observe that maxima of the density 
coincide with maxima of g{x) and minima of V{x), cf. Figs. [6]and[7l 

A stability region, summarizing the results of many runs of the simulations, can also be identified in this case, as 
shown in Fig. [T^ cf. the stability diagram shown in Fig. 2] for the cn solutions (without the linear potential). The 
dashed curve with arrows designates the existence border for the solutions of this type, as determined by condition 
> 0. It follows from Eq. p6)) . that, in the case shown in Fig. [121 this existence condition amounts to fc^ < 



VI. SOLITONS 



A. Stationary solutions for bright solitons 



In the limit of A; —> 1, all the above types of the exact periodic solutions go over into localized structures (solitons) 
pinned by the corresponding localized pseudopotential, possibly in the combination with the localized (trapping) 
linear potential. In this limit, the solutions of both the cn and dn types, as given by Eqs. ([5]), ([5]), dH]), or ([27]), 
([30)l . reduce to a common expression for the bright soliton. 



.go(l + 2b) 

Ao\{x) = W — -2 — — , (50) 
y cosh X + b 

with fiso\ ~ ^1/2- The nonlinearity-modulation profile which supports this solution, in the absence of the linear 
potential, is given by Eqs. Q and ([7]), in which cn{x) is replaced by sech x, i.e., 

, ^ go+gisech^x 
1 + sech (x) 

and coefficient gi is replaced by its limit form for k — I, i.e., as follows form Eq. ([7]), 

g,^k^l)^'4^. (52) 

^ 2 26+ 1 ^ ^ 
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FIG. 11: The same as in Figs. [6l[71 and[9l but for a typical stable solution of the sn type found in the presence of the linear 
potential, for go — —1, b — —0.8, r = 1, gi — 2, and k = 0.7. 
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FIG. 12: The region of the stability for the sn-type periodic stationary solutions, supported by the linear potential, as found, 
for go — —1 and gi = 2, in the plane of parameters b and k, which are defined in Eq. (|35|l . cf. Fig. 4. 



The soliton may be stable only if g(x) given by Eq. (|5ip has a maximum at a; = [for a minimum of g(x) at 
a; = 0, an obvious variational consideration of energy (j3|) immediately predicts an instability of the soliton placed at 
the maximum of the respective pseudopotential] . The general condition for the existence of the maximum of g{x) 
at a; = was actually obtained above, gi > ggb. In the present case, with gi taken as per Eq. ([5^ . it amounts to 
— 1 < 6 < 0, which is, thus, a necessary condition for the stability of the soliton. Since solution (|50p exists only for 
go {1 + 2b) > 0, we finally conclude that the solitons may be stable only in the case of 



go^-l, -l/2<6< -1. 



(53) 
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In fact, the nonlinearity is sign-changing in this case, because g{x = ±00) = go = —1, while g{x = 0) = (1/2) (2 — 
|5|)/(2|5| - 1) > 0. Lastly, the norm of soliton §UI is 

A more general bright-soliton solution can be obtained by setting fc = 1 in the cn or dn solutions including the 
linear potential, i.e., Eqs. Q, (O, (ITSt . (IT6l) . (IT7|) . and (ITSl) [recall that r cannot be set equal to 1 in this case, and 
gi is a free parameter, while Eq. ([7]) is irrelevant]. This leads to /i = 1 — r^/2 and 



y 2 {go- gi/b) [cosh (ra;) + b\ 



The existence condition for solution ((55|) is 

50 > 51/^, (56) 

which is different from that at which soliton (jSOp exists. Further, the nonlinearity-modulation function supporting 
solution (I55p is given by the same expression (jSip as above, with the difference that gi is a free parameter in it, and 
the respective form of the linear potential is 

1 + sech [rx) 



Vo = b+- 



l-5 + 3ii±^ 

9ob - gi 



(58) 



Potential ([57]) stabilizes the soliton if it represents a potential well at x = 0, which means Va < b. A straightforward 
analysis of expression ([55)1 . taking into regard Eq. ([5^ . demonstrates that the latter condition is met for 

go(6-l) > (2/fe)(26+l)gi. (59) 



B. The stability of bright solitons 

Systematic simulations of Eq. ([T|) with V{x) — 0, starting from the initial condition corresponding to a perturbed 
bright soliton, demonstrate that the solitons are indeed stable in region ([5^ . where this should be expected, as argued 
above. An example of the stable soliton is shown in Fig. [131 For values of b approaching —1, see Eq. (|53| . the initially 
imposed random perturbations gradually grow, which, however, is interpreted as the growing sensitivity of the solution 
to the random perturbations, rather than as a true instability. A technical definition of the "low sensitivity" to the 
random noise may be adopted in the form of the requirement that the relative amplitude of the perturbations must 
remain below 3%, after the evolution time t = 5000, if the initial perturbation amplitude was 1%, as fixed above. In 
particular, for the same parameters as in Fig. 1131 r — 1 and gi = 2.4, the so defined sensitivity border is found at 
b — —0.75. For instance, the final perturbation amplitude is 2.82% and 3.36%, at 6 = —0.7 and —0.8, respectively. 

In the presence of the linear potential, stable solitons can be found both for go = —1, as above, and for go — +1, 
see Figs. [MlfTHl For instance, fixing go — +1, r = 1 and gi = —2, as in Figs. [T4l and [T5l we find that the solitons 
are stable for < 6, and unstable for 6 > 6 [the soliton is unstable at 6 = 6, as shown in Fig. [121 but still stable-at 
least, up to t = 5000 -at b = 5.9 (not shown here)]. In fact, the transition from the stability to instability, in the form 
of a spontaneous escape of the trapped soliton, as observed in Fig. [151 is explained by the competition between the 
self-repulsive localized nonlinear pseudopotential, and the attractive linear trapping potential, see panels (b) and (c) 
in Figs. [HlandfTSl Note that parameters corresponding to both figures satisfy condition (jS^]) . under which the linear 
potential is the trapping one. 

Lastly, the combination of go = ~1 with the linear potential makes both the nonlinear pseudopotential and linear 
potential attractive, hence the solitons are definitely stable in this case, see an example in Fig. [THl For the parameters 
chosen as in this figure, go = —1, r = 1, 171 = 2, the above-mentioned technical definition of the sensitivity of the 
trapped soliton to the random noise (the perturbation amplitude at the level of 3% by t = 5000) sets the "sensitivity 
border" at 6 = —0.65, the solitons being "hypersensitive" at —0.65 < 6 < — 1. 
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FIG. 13: The same as in Fig. (2] but for a case of a stable bright soUton, obtained as a limit of the cn solutions with k — 1 
[accordingly, panel (b) displays the localized sign- changing nonlinearity-modulation function g(x) which supports the soliton]. 
Other parameters are go = —1, b — —0.6, r — 1, gi = 2.4. 




C. Dark solitons 



The limit form of the exact solutions of the sn type for fc = 1 are dark solitons. In the general case, when the 
respective linear potential (|44)) is included, the dark-soliton solution can be obtained from Eqs. ((36)) and p6|) in the 
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FIG. 15: The same as in Fig. 1141 but for 6 = 6. In this case, the nonhnear self-repulsive pseudopotential turns out to be 
stronger than the linear trapping potential, hence the bright soliton is unstable, being spontaneously expelled from the bound 
state. 



following form: 



'0dark-sol(a;) 



36 



r(l + b) tanh (rx) 



(60) 



The nonlinearity-modulation function and linear potential which support this dark soliton are obtained from Eqs. 
([55)1 and (gH) by setting k = I in them: 



(Uo) 



9{x) 
V{x) 

dark— sol 



(.90 + gi) - .gisech^ {rx) 
{l + h)-b sec\\ {rx) 
(^o)dark-soitanh^ {rx) 
{l + h)-h sech^ {rx) 



1 - 6 (2 + 3&) + 



3gi (1 + b)' 
{gob - gi) 



(61) 
(62) 
(63) 



The stability of the dark solitons will be considered elsewhere. 



VII. CONCLUSIONS 



This work reports three large classes of exact periodic solutions for the one-dimensional GPE (Gross-Pitaevskii 
equation) with the pseudopotential represented by the periodic modulation of the nonlinearity coefficient, which can 
be created in the BEG by means of the Feshbach-resonance technique. In the general case, the model includes a 
periodic linear-lattice potential too. The modulation function of the corresponding NL (nonlinear lattice), the linear 



19 




X X 

(b) (c) 

FIG. 16: The same as in Fig. [14] (also a stable soliton), but for go — —1, b = —0.5, r = 1, gi = 2. Panels (b) and (c) 
demonstrate that both the nonlinear pseudopotential and linear potential are attractive in this case. 

potential (if any), and the exact solutions are expressed in terms of the Jacobi's elliptic functions of three types-cn, 
dn, and sn, which give rise to the three different classes of the solutions. In the absence of the linear potential, the 
solutions depend on two continuous free parameters, b and fc, and the sign parameter, gg. If the linear potential is 
included, full solution families feature two additional continuous parameters, gi and r. The stability of the exact 
solutions was tested by means of systematic direct simulations of the perturbed evolution, except for dark solitons, 
whose stability will be considered separately. 

Subfamilies of the pseudopotentials and respective exact solutions have been identified with both sign-changing 
and sign-constant nonlinearity-modulation functions. Density maxima of the solutions may coincide with minima or 
maxima of the periodic pseudopotential; in the former case, the exact solutions may represent the ground state of the 
BEC. In the absence of the linear potential, only the solutions with the density maxima collocated with minima of 
the pseudopotential may be stable (although this condition alone is not sufficient for the stability, in most cases). On 
the other hand, the addition of the trapping linear potential may stabilize solutions whose density maxima coincide 
with maxima of the nonlinear pseudopotential. 

In the limit of the infinite period, the exact solutions reduce to solitons. The solitons are stable in a large part 
of parameter regions where they exist. In particular, they may be stable if located at a maximum of the nonlinear 
pseudopotential, provided that the competing trapping linear potential is strong enough. 

This work calls for continuation in other directions. In the models considered here (with both the periodic and 
localized pseudopotentials), it is interesting to find, by means of numerical methods, a full set of nonlinear states, 
which should reveal how the exact solutions are embedded into the full family, and exactly identify the respective 
ground state. Further, it is desirable to investigate the dynamical stability of the solutions in a rigorous form, through 
the computation of eigenvalues for modes of small perturbations. Another interesting issue is a possibility of finding 
similar families of exact solutions for systems of two or three coupled GPEs describing multi-component BECs. In 
fact, the methods used in earlier works [l^, [l^,!!^] for obtaining families of exact solutions in coupled NLSEs can be 
used in the latter context. 
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